Секвенирование транскриптома

1. Анализ качества чтений

Команда Описание
fastqc chr2.1.fastq Выдает характеристики качества сырого чтения
java -jar /nfs/srv/databases/ngs/suvorova/trimmomatic/trimmomatic-0.30.jar SE -phred33 chr2.1.fastq chr2.1_trim.fastq TRAILING:20 MINLEN:50 Команда TRAILING удаляет буквы с качеством ниже установленного с конца рида. Команда MINLEN удаляет риды длиной меньше установленного значения
fastqc chr2.1_trim.fastq Выдает характеристики качества чтения

Сверху представлены качества чтений соответственно до и после триммированя, как видно качество чтений и так было хорошим поэтому оно особо не изменилось

2-3. Картирование чтений. Анализ выравнивания

Команда Описание
hisat2-build chr2.fasta chr2 Команда hisat2-build индексирует референсную последовательность ДНК
hisat2 -x chr2 -U chr2.1_trim.fastq -S chr2.1.sam --no-softclip Команда hisat2 выравнивает прочтения и референс. -x указывает на проиндексированный референс, -U - на файл с ридами, --no-softclip запрещает убирать нуклеотиды с концов ридов. Параметр --no-spliced-alignment запрещает разделять риды и картировать их на разные места генома, а так как в данном случае возможно картирование разделенных ридов на разные места генома в результате сплайсинга, этот параметр указывать не нужно.
samtools view chr2.1.sam -b >> chr2.1.bam Переводит файл с выравниваниями в бинарный формат
samtools sort chr2.1.bam sort Samtools sort сортирует полученный бинарный файл по начальным координатам
samtools index sort.bam Samtools index индексирует файл
Выровнялось 12345 рида, не выровнялось 54, выровнялось по 2 раза 0 ридов

4-5. Подсчет чтений. Анализ результатов

Команда: htseq-count -f bam -s no -m intersection-nonempty -i gene_id sort.bam /P/y14/term3/block4/SNP/rnaseq_reads/gencode.v19.chr_patch_hapl_scaff.annotation.gtf >> count2.txt

-f отвечает за формат входного файла с выравниванием (.sam или .bam; .sam по умолчанию)

-s - цепь, на которой ищутся гены с ридами. Ставим no, так как, анализируя транскриптом, мы не знаем, с какой цепи получается транскрипт

-i - какой аттрибут использован для характеристики места картирования рида (варианты gene_id/transcript_id)

-m - для задания перекрытия ридов, берем intersection-nonempty, тогда риду присваивается gene_id в случае хоть какого-то непустого перекрывания

Выдача скрипта:

ENSG00000115053.11 12058

ENSG00000202400.1 5

ENSG00000206885.1 12

ENSG00000233538.1 2

__no_feature 268

__not_aligned 54

12340 чтениq не легли в границы генов. ENSG00000115053.11 - ген нуклеолина, ядрышковый белок принимающий участие в создании рибосом. ENSG00000233538.1 - ген неописанного белка